num_cores = str2num(getenv('SLURM_CPUS_PER_TASK'));
c = parcluster;
poolobj = parpool(c, num_cores);

clear all
clc

cd('../cvx_2.2.2')
cvx_setup

warning('off','MATLAB:nargchk:deprecated')
warning('off','MATLAB:illConditionedMatrix')

cd('../experiment')

data = readtable('data.csv');

J = 3;    % # of goods

data_wide.choice = data.choice;
data_wide.x = [data.x1 data.x2 data.x3];
data_wide.z = [data.z1 data.z2 data.z3];

data_wide = struct2table(data_wide);

[data_trans] = transform_data_01_exp(data_wide,J);

m_regr = [3 4 3*ones(1,J-2) 4 4 3*ones(1,J-2)];

[combos_all_unique] = combinations_10(m_regr,J);

% constraints

[Aineq,bineq] = constraint_mat_theta_2(combos_all_unique,J);

Aeq = [];
beq = [];

ntheta = size(Aineq,2);

N = size(data_wide,1);
[theta_NPD,fval_NPD,flag_NPD] = estimation_NPD_cvx(data_trans,J,N,m_regr,Aineq,Aeq,bineq,beq,...
        combos_all_unique,ntheta);

[ratio_trimmed_mean,mean_ratio,trimmed_mean_ratio] = script_ratio_rand_exp_fun_average(theta_NPD,m_regr,data_wide,J,ntheta,combos_all_unique);
ratio_trimmed_mean

clear data* A* b* der_den der_num
    
save('estimates.mat')
